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ABSTRACT 

The significance of hierarchical clustering on the density profile and mass- 
temperature scaling relation for galaxy clusters is examined using hydrodynamic N- 
body simulations. Clusters formed hierarchically are compared with clusters formed 
with the initial density fluctuations on sub-cluster scales removed via smoothing. 

The universal profile, as described by Navarro, Frenk, & White, is not a by-product 
of hierarchical clustering. It is found to fit the mean profiles of clusters formed both 
hierarchically and otherwise. The Hernquist profile is also found to fit the data well. 
The characteristic radius, r^, moves outward from 0.1ii200 to 0.2i?200 when the initial 
substructure is eliminated. Interior to r^, pdm oc r~^'®, regardless of initial smoothing. 
Exterior to this radius, the profile is marginally shallower in the non-hierarchical case, 
with pdm oc r~^'^ compared with pz^A/ oc r"^'''. 

The mass-temperature scaling relation maintains the form T (x M regardless of 
cluster formation method. The normalisation varies at the 20% level, which is at the 
level of the intrinsic scatter, with the non-hierarchical simulations producing the cooler 
clusters. 

Subject headings: hydrodynamics - methods: numerical - galaxies: clusters: general - 
large-scale structure of Universe - X-rays: galaxies 



1. INTRODUCTION 

The largest virialised objects in the universe are galaxy clusters. On scales larger than this, 
the quasi- linear collapse of large scale structure dominates. Tthese clusters were initially identified 
optically from maps of galaxies which indicate galaxies tend to be found in together in groups, some 
of which contain many hundreds of galaxies. Subsequently, these clusters of galaxies were found to 
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contain a diffuse intergalactic medium composed of hot (10^ K « 10 keV ) gas. This gas is luminous 
{rsj 10^^ L0 ) in the x-ray regime and, consequently, observable from above the earth's atmosphere. 
Two such satellites in particular, Einstein and ROSAT, have been sources for catalogues of x-ray 
luminosities (Gioia et al. 1990) as well as luminosity maps with sufficient resolution to map the 
coarse structure of these clusters (McMillan, Kowalski, &; Ulmer 1989; Breen et al. 1994). This hot 
gas constitutes the bulk of the baryonic matter in galaxy clusters. 

Satellites with instruments able to take spectral data in the x-ray band, most notably the ASCA 
satellite, have coarsely resolved the spatial temperature distribution of this gas (see Markevitch & 
Vikhlinin (1997), for references). A new breed of x-ray telescopes including Chandra, XMM, and 
Astro-e will be able to better map the temperature and mass distributions of clusters (Schindler 
1999), potentially to a level which can discriminate competing cosmological models. 

Being so luminous, galaxy clusters are observable to great distances and, consequently, at ear- 
lier stages of their evolution. This permits tests of many cosmological parameters (Henry et al. 
1997; Donahue et al. 1998, 1999; Gioia et al. 1999). Galaxy clusters are observable in ways other 
than their x-ray signature and population of galaxies. In the microwave spectrum, there is the tem- 
perature decrement imposed on the cosmic microwave background via the Sunyaev-Zeldovich effect 
by which electrons in the hot gas alter the energy spectrum of the CMB photons (see Carlstrom, 
Joy, & Grego 1996; Komatsu et al. 1999; Andreani et al. 1999, for examples). Optical imaging 
of the clusters also provides independent measurements of the total mass of the clusters via the 
gravitational lensing of background galaxies (see Clowe et al. 1998; Makino &: Asano 1999, for 
examples) , though the method still has details to be clarified (Metzler et al. 1999) . 

In the standard model of cosmological structure formation, clusters are the most recently 
relaxed large-scale objects. Consequently, though they are relaxed, they potentially carry with them 
information related to their formation. Larger-scale structure is in a state of linear or quasi-linear 
evolution and, hence, still responding to density perturbations that have existed since early times. 
Smaller virialised objects such as galaxies and stars individually contain little or no manifestation 
of the character of the initial density perturbations from which they formed. Having had ample 
time to dynamically relax, their present evolution is dominated by the physics of stellar evolution. 
Clusters are representative of structure on the intermediate scale. 

There is ample evidence that the dominant form of matter in the universe is unseen and not 
baryonic in nature. This dark matter is fundamental to the standard model. Though the presence of 
dark matter may be inferred from observational data, by its nature it cannot be presently readily 
observed directly. Numerical simulations which model this coUisionless dark matter provide a 
powerful tool to achieve the goal of understanding the history and dynamics of structure formation. 
Two recent results stemming from numerical simulations have added to our understanding of the 
morphology of galaxy clusters as well as provided a tool for inferring the total mass of clusters 
from their measured temperatures. These are the universal profile of dark matter and the mass- 
temperature scaling relation. 
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1.1. A universal density profile 

The works of Navarro et al. (1995) and Navarro, Frenk, &: White (1996) have built a case for 
a common density profile for the dark matter in galaxy clusters. This density profile is found to fit 
clusters spanning a large range of masses. It is contended that this 'universal' profile has only one 
free parameter corresponding to the density at cluster formation (Navarro, Frenk, &; White 1997). 
However, the form of the profile is derived from numerical simulations of hierarchically formed 
clusters. There is evidence that it is the hierarchical nature of the cluster formation itself that is 
responsible for the universal profile (Syer & White 1998). Determining the validity of this would 
be valuable to understanding the true breadth of its universality. 

1.2. The mass-temperature scaling relation 

Estimations of the total mass of the clusters have been made using the observed x-ray luminos- 
ity (and hence, baryon) distribution and the assumption of hydrostatic equilibrium using so-called 

/3- model fits (Fabricant & Gorenstein 1983; Jones &: Forman 1984). Evrard, Metzler, &; Navarro 
(1996) describes a correlation between the total mass of clusters of galaxies and the mass-weighted 
mean temperature of the hot, x-ray emitting gas in the interior of the cluster. This would provide 
an independent method of measuring the total masses of the clusters, requiring only a measurement 
of the temperature of the gas. This relationship is a natural result of the hydrostatic state of the 
gas. However, derivation of this relationship uncovers a dependency on both the gas density and 
temperature profiles. Hence, a sensitivity of these profiles to the cosmogony would weaken the 
relation's utility. Understanding the behaviour of the mass-temperature scaling relation on the 
cosmogony is essential to the confidence level put in the application of the relation to real data. 

1.3. Hierarchical clustering 

In hierarchical clustering, the largest structures forming at a given time do so via the amalga- 
mation of many smaller structures which have formed at an earlier time. This is owing to the form 
of the initial density perturbation spectrum in which small scale perturbations have higher initial 
amplitudes than large scale. In contrast, non-hierarchical clustering involves structure formation 
from the collapse of large structures with smooth density distributions. Though the details are 
still not clear, the results of numerical simulations compared with observations support the theory 
that we live in a universe in which structure is formed hierarchically. The degree to which the 
hierarchical nature affects galaxy clusters is not entirely clear. 

In this paper, the issue of whether the mass-temperature scaling relationship and the common 
density profile are universal between hierarchical and non-hierarchical cluster formation scenarios 
is explored. By using these extreme cases, the significance of hierarchical clustering itself will be 



-4- 



determined. Through inference, their dependence on the presumed cosmological model will also 
be illuminated; if they are found to be independent on hierarchical clustering, then it would be 
expected that the precise form of hierarchical clustering would also be irrelevant. 

The analysis will concentrate not on the properties of individual clusters, but on the global 
mean properties of scaled quantities. 

The layout of the paper is as follows. In Sec. 2, the simulations arc described. A brief 
description of the analysis methods is given in Sec. 3. Verification that the clusters extracted 
arc in hydrostatic equilibrium is made in Sec. 4. The dark matter density profile is examined in 
Sec. 5. The derivation of the mass-temperature relationship and the results relevant to clusters are 
presented in Sec. 6. The conclusions are summarised in Sec. 7. 
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2. SIMULATIONS 
2.1. Numerics 

The N-body AP^M- SPH code, HYDRA (Couchman, Thomas, & Pcarcc 1995) was used for all 
simulations. This is a multi-level Adaptive Particle-Partiele Particle-Mesh (AP'^M) N-body code 
with gas dynamics simulated by the Lagrangian Smooth Particle Hydrodynamics (SPH) method. 
The gravity is calculated using a particle-mesh scheme for the large scale gravitational fields. Short- 
range forces are calculated by summing particle-particle forces. Isolated regions of high number 
density have sub-grids adaptively placed around them, allowing the efficient particle-mesh method 
to be used locally. As such, it is well suited for the large dynamic range of scales involved in cluster 
studies. A comparison of the leading hydrodynamic codes, including HYDRA, designed for studying 
cosmological scenarios, finds general consistency among the codes (Frenk, White, Bode, Bond, 
Bryan, Cen, Couchman, Evrard, Gnedin, Jenkins, Khokhlov, Klypin, Navarro, Norman, Ostriker, 
Owen, Pearce, Pen, Steinmetz, Thomas, Villumsen, Wadsley, Warren, Xu, &; Yepes 1999). 



2.2. Cosmology 

All simulations assumed the flat cosmology given in Table 1. The matter was evolved in a 
box with 40/i~^ Mpc sides in co-moving coordinates, with the Hubble expansion constant. Ho = 
/ilOOkm/s . This permits the formation of a sufficient number of clusters for statistical purposes. 
A box of larger dimensions would impose too great a penalty on the spatial and mass resolution of 
the simulation. A smaller box would lack the long-wavelength perturbations required to produce 
large clusters. The simulation had periodic boundary conditions. 

Cooling is neglected. The cooling time for the bulk of the cluster gas is estimated to be well 
over the age of the universe. The peak luminosity is on the order of 10^^ erg s~^. The thermal 
energy found in the gas of a system containing 10^^ M0 of gas at 5 keV is on the order of 10^^ erg 
giving a cooling time, tcool = 10^^ s 3 x 10^^ yrs . There is insufficient resolution due to the 
limit set by the gravitational softening parameter to properly model the cooling flows which are 
inferred in the inner 200 kpc of clusters (Allen & Fabian 1997). It would require a great deal more 
resolution in order to model these flows correctly owing to the overcooling induced by the inability 
of SPH to properly calculate the gas densities in the cores of clusters. As well, when modelled with 
sufficient resolution, the cooling flows simulated are much greater than those observed, likely due to 
a lack of feedback mechanisms such as supernova energy input which reheats the gas (Suginohara 
&; Ostriker 1998). In any case, the amount of gas inferred to have cooled in these flows is on the 
order of 10^^ Mq , which is not significant to the larger scale distribution of gas in the cluster. 
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2.3. Initial conditions 

The initial density perturbations were established by displacing the particle positions from 
a uniform cubic grid using, in the standard way (Efstathiou, Davis, Prenk, & White 1985), the 
Zel'dovich approximation for growth of density perturbations in the linear regime (Zel'dovich 1970). 
This method first creates a representation of the density perturbations in Fourier k-space using a 
supplied power spectrum form supplemented with Gaussian random fluctuations. This is then 
transformed into an initial density field. The gravitational forces felt by the set of particles, dis- 
tributed at the nodes of the uniform mesh, by this density field are calculated. The particles are 
then displaced in the direction of their respective forces an amount proportional to the force. This 
produces a distribution of particles with a density field following the previously calculated density 
field. 

The initial redshift for the simulations is ZinUial = 75. This was chosen to keep the maximum 
displacement incurred during the establishment of the density field to less than V2 the initial grid 
spacing. This makes negligible the errors incurred in using the linear Zel'dovich approximation to 
a non-linear system. The initial power spectrum of the density fluctuations follows a power-law of 
n = -1, P{k) oc k~^. 

In order to ascertain the effects hierarchical clustering has on the matter, runs with initial high- 
spatial-frequency density perturbations were compared to runs in which these density perturbations 
were suppressed. The perturbation suppression was accomplished via two methods. 

For the first, the initial density field was convolved with a tophat. This modification to 
the power spectrum takes the form P'{k) = P{k)W{kr smooth) where W{x) is the tophat window 
function given by, 

, 3(sina; — xcosx) 
= ^3 • (1) 

A choice of r smooth = 7h~^ Mpc for the smoothing length was selected to suppress those perturba- 
tions of spatial size less than the size of a cluster formation region. Consequently, the clusters in 
this scenario form from the uniform collapse of structures of a size on the scale of or larger than 
those that formed the clusters seen at the present. To span the regimes, a set of initial conditions 
smoothed with r smooth = 3/i~"^ Mpc was also evolved. In this case, the early generations of struc- 
ture formation would be suppressed, but the largest objects at the final epoch would have formed 
from the merger of a few objects. 

The second method for smoothing the initial power spectrum is a low-bandpass filter applied 
to the power spectrum of the density-perturbations before the particles are displaced. This filter 
is of the form 

P'{k) = P{k)e y ) , (2) 



where n is some positive integer which controls the sharpness of the cutoff which was set to n = 16. 
Two sets of initial conditions for which r smooth = '^h~^ Mpc and 14h~^ Mpc were created using 
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this method. The power spectra for the initial conditions created with the top-hat smoothing as 
well as the frequency cutoff are illustrated in Fig. 1. 

The first minimum in W{x) is at A; | ^ ^^t^ - Hence, tophat smoothing over a radius 
of 7h~^M.Tpc removes more power at frequencies k < ^ir/rsmooth than a frequency-cutoff filter 
using the same value of r smooth- The significance of this is made apparent by noting that the 
frequency-cutoff filter is essentially a convolution in real space of the density distribution with the 
window function, W{x). The W{x) function can crudely be thought of as a top-hat with a radius 
given by the position of the first minimum. This implies a frequency cutoff of 2'K{7h~^ Mpc)~^ is 
roughly comparable to smoothing with a top-hat of radius bh~^ Mpc . Our choice of r smooth for 
the frequency filtering reflects this, with a frequency cutoff of 27r(7/i~^ Mpc)~^ producing an initial 
density distribution roughly matching that produced by "ih^^ Mpc top-hat smoothing. Similarly, 
the cutoff of 27r(14/i~^ Mpc)~^ produces results similar to the 7h~^ Mpc top-hat. This is verified 
at a visual level in the distribution of matter at the end of the simulations (Fig. 2). For brevity, 
subsequent references to the scale of frequency cutoffs will omit the 27r term, which may be assumed. 

Since the top-hat smoothing removes power at all frequencies except 0, the power on all scales 
was increased during the 3^~^ Mpc top-hat smoothing to maintain the same RMS fluctuations on 
scales of 8/i~^ Mpc (as) as found in the unsmoothed set of initial conditions. This was not done for 
the (14/i~^ Mpc)~^ frequency filtering, as this removes power up to scales larger than Sh~^ Mpc, 
nor was it done for the 7h^^ Mpc top-hat filter for which first minimum in k-space is too close to 
the scale of {^h~^ Mpc)^^. Frequency filtering over (7h~^ Mpc)^"*^ does not require this correction. 

To examine the effects of resolution, the runs were performed at two resolutions. In one set, 64^ 
particles each of gas and of a collisionless dark matter component were evolved. In the other set, 
32^ particles of each were used. The gravitational softening was 0.02/t~^ Mpc and 0.04^~^ Mpc for 
the respective runs. The details of the simulations are given in Table 2. The difference in structure 
at t = 1 (i.e., the present) is clearly illustrated in Fig. 2 which maps the projected gas density in 
the each of the simulated volumes. 



2.4. Two-body heating 

The effect of two-body interactions is minimised by the gravitational softening parameter, e. 
Without a gravitational softening parameter, gravitational interactions heat the gas via an exchange 
of energy between the dark matter and gas phases during short-range encounters of pairs of these 
particles (Steinmetz k. White 1997). If the mean interparticle spacing is on the order of or greater 
than e, then the time, t2-body, required to heat the gas to its temperature, T, via this process is 
given by, 

iw^) = — } — ^f7^V/3- - 

\1U yrs / /lnA\ / moM \ mDM \ ' 



'\nA \ ( mpM \ ( mpM \ ( Pdm \ U2 

K 5 ) \^1O1»M0 ) \ rUgas ) \ Po J 
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where In A is the Coulomb logarithm (3-7 for most simulations), uidm and mgas are the respective 
masses per particle for the dark matter and gas, pdm is the dark matter mass density, and pc is 
critical density. 

Without a softening term in the gravitational forces, two-body interactions would heat the gas 
in the centre of simulated clusters to 10^ K in 0.5 Ga . Indeed, it can be shown for the simulations 
presented here, in the absence of a gravitational softening term, two-body heating would be sig- 
nificant for radii less than ~ O.lTMpc for the simulations using 2 x 64^ particles and ~ 0.5 Mpc 
for those using 2 x 32^. In the dense cores, where this process is most significant owing to the 
dependency, t2-body oc pj)\j, there are on the order of 70 particle within a softening length, e, of 
any other particle. This obliviates the problem, there. As the density drops, however, this number 
drops. Again, using the density profile information described in Sec. 5.2, it may be found that the 
mean number of particles within e falls to unity at a radius of ~ 0.1 Mpc , independent of the spatial 
resolution of the simulations, since e scales with spatial resolution in these simulations. This implies 
that two-body heating is probably a factor in the regimes of 0.1 to O.lTMpc for the high-resolution 
simulations and 0.1 to 0.5 Mpc for the lower. In the outer halo, poM falls off faster with radius 
than does the factor, T^^^, so two-body heating does not become a factor except, perhaps, for the 
cold void particles. 
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3. ANALYSIS 

3.1. Cluster selection 

Clusters were selected by a tomographic deprojection method. The gas distribution, extended 
over its smoothing radius, was projected on to three planes: the x-y, the x-z, and the y-z. The 
resolution of these projections is ~ 100 kpc . Peaks were found in these three planes which were 
then tomographically deprojected to recover the (x,y,z) coordinates of the density maximums. The 
deprojection method involves the following steps: 

1. Match the peak positions in the x domain (within some tolerance) for the x-y and x-z pro- 
jections to find a list of possible y-z positions. 

2. Match the y-z positions postulated from the previous step with the y-z positions found from 
the peaks in the y-z projection. 

3. For each of the matches, use the x position from step 1) to get an x-y-z position. 

This method is comparable to projecting the gas densities onto a 3-dimensional mesh of resolution 
(box size)/L, but requires 3L^ elements of information instead of Li^ which, beyond the memory 
requirement savings, speeds up the peak-finding algorithm. 

Cluster searching using the friends-of- friends method (fof vl.l ^) did not give satisfactory 
results. This method links together particles closer than a distance given by a linking length and 
then associates to a particle all its neighbour particles as well as their neighbours and so on. The 
association is then called a cluster. Close clusters are often connected by bridges allowing the 
association to 'percolate' through. As such, the grouping of cluster members was found to be 
too sensitive to the linking length parameter. For a choice of a linking length suitable for the 
largest clusters, the smallest clusters were not properly recovered. Another popular routine, SKID ^ 
(based on denmax), was tested. This routine 'freezes' the particles, then allows them to gravitate 
together to form tighter groups, essentially severing the bridges between clusters. A friends-of- 
friends procedure is then allowed to create the list of associations. Though the routine denmax 
performed better than FOF, particularly for the 2 x 32^ simulation, the amount of CPU time required 
for the high-resolution simulation was unacceptable, taking days instead of hours. 

The centres of the density peaks are further refined by cutting out a small sphere of particles, 
centred on the estimated density peak, and finding the mean position of the particles with the 
highest densities calculated previously by the SPH algorithm. 



^see http: / /www-hpcc. astro.washington.edu / tools/FOF/ 
^see http://www-hpcc.astro.washington.edu/tools/SKID/ 
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Overdensity radii, Rg, were calculated for the clusters for overdensities of ^ = ^ = 200 and 
500. An overdensity radius is simply the radius from the centre of a cluster within which the ratio 
of the mean density to the critical density of the universe is equal to some value, referred to as the 
overdensity. Though spherical symmetry is not required, the derivation holds more relevance for 
spherically symmetric systems. The expected radius of virialisation derived from analytic models 
involving simple spherically symmetric collapse corresponds to that of an overdensity of 178 (i.e., 
i?i78). However, the overdensity radius, Rg, is not particularly sensitive to 5. The overdensity 
within the calculated overdensity radius is accurate to a factor of where N is the number of 
particles within Rg. With the overdensity radii, there is a corresponding overdensity mass, Mj, 
which is the mass contained within 

A lower limit to the size of the clusters was set by the requirement that each cluster, within 
the overdensity radius of -R50O) have at least 300 gas particles and 300 dark matter particles. This 
ensures the densities are calculated correctly (see Appendix A). 

The details of the sets of clusters found using this method are given in Table 3. 

A similar procedure was done using the dark matter mass distribution. The tomographic 
dcprojcction step for finding the clusters first required SPH-likc density estimates be made for the 
dark matter and appropriate smoothing lengths be calculated in order to project the masses onto 
the x-y, x-z, and y-z planes. Little difference was found in the positions or number of the final 
clusters selected when the dark matter distribution was used instead of that of the gas. 

3.2. Cluster profiles 

Profiles of various parameters were calculated for the clusters by summing the contributions of 

the particles falling in radial bins centred on the clusters. The bins were separated exponentially. 
Found for each of the dark matter and gas components were the number of particles, the volume- 
weighted mean density, and the mass-weighted mean thermal energy. The innermost bin was 
selected to be at a radius greater than that which encloses 300 particles to ensure that the gas 
dynamics well-models the density gradient (see Appendix A). 

The thermal noise of the SPH particles provides a negligible contribution (1%) towards the 
thermal energy of gas particles in the haloes. This contribution is ignored in all cases. 

Mean profiles were calculated using the radially binned profiles described previously. The 
individual cluster profiles were scaled radially by the overdensity radius, R200, as well as weighted 
by their respective cluster mass. For the density profiles, this weighting was found to have only a 
small effect, however, validating the use of mean profiles. 

Density and temperature profile forms were fit to either the mean profiles or to the particle 
distributions of density or temperature versus radius. Details of the fits are given in Sections 5.2, 6.4 
and 6.5. In all cases, the fits were found by minimising the with the variance derived from the 
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distribution of values in radial bins. Again, the distances were scaled by the overdensity radius, 
-^200- When fitting to the particle distributions of density or temperature versus radius, all the 
data from selected clusters in a given simulation were used simultaneously. 

Errors for the fits to the density or temperature versus radius distributions were determined 
using a bootstrap method. In this method, one half of the data arc fit at a time, giving a set of 
parameters to the fit. This is done a multitude of instances for the same data, selecting a separate 
random sample of points in the distribution for each instance. These fits produce a set of values 
for each parameter to the fit, from which a mean and deviation is derived. 

It is the SPH-calculated density that is used when the particle distribution of density versus 
radius is fit. Both the dark matter and gas particles had their densities calculated more precisely 
by performing an SPH summation over a value for the smoothing parameter, h, such that 2h 
encompasses exactly NgpH = 32 particles. For the gas, this primarily affects the low and high 
density regimes since during the simulation there is a lower and upper bound set on h by the 
gravitational softening length and the spacing of the mesh used to facilitate neighbour searching, 
respectively. The densities were calculated without a self-density term to provide a volume- weighted 
average local density at each of the particle positions. 
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4. HYDROSTATIC EQUILIBRIUM 



The derivation of the mass-temperature scahng relation, as weh as the /3-model method for 
estimating the mass of galaxy clusters require that the clusters be in hydrostatic equilibrium within 
some distance from the cluster centre. This can be tested in the numerical simulations. Integration 
of the equation of hydrostatic equilibrium produces an estimate for the gas pressure which can be 
compared with the actual gas pressure profile. This integral is given by, 



where M(r) is the total mass interior to the shell of radius r and {pg{r)) is the mean gas density at 

r. The calculated value, Pg{r), is the gas density at r required to satisfy the condition of hydrostatic 
equilibrium. In practice, integration to r ^ oo is both impractical, due to the finite box size, and 
inappropriate, since the cluster is not isolated. Instead, a radius at which the presence of nearby 
clusters influences the profiles is chosen. This is always at a distance such that the contribution to 
the integral will be minimal. 

The analysis described above was performed on the sample of cosmological clusters described 
in Sec. 2. The results for the hierarchically formed clusters, shown in Fig. 3, indicate that the 
clusters have mean radial pressure profiles that are consistent with the condition of hydrostatic 
equilibrium. This is particularly true within the virial radius indicated by i?200- The same is true 
for clusters formed non-hierarchically (Fig. 4), for which we used the data for the clusters formed 
from initial conditions smoothed by a 7h^^ Mpc top-hat as the sample. Consequently, it may be 
accepted that the clusters in this study are approximately hydrostatic out to their virial radii. 




(4) 
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5. UNIVERSAL PROFILE 

Numerical simulations have indicated that, for a Cold Dark Matter (CDM) model, the dark 
matter density profile of galaxy clusters may be fit by a form with only one free parameter set 
by the mean density of the universe at the time of collapse (Navarro et al. 1996). The existence 
of STich a universal density profile (or NFW profile) would be a useful tool as well as a test for 
theories of cluster formation. The applicability of the universal profile to clusters formed in other 
cosmological models would strengthen its use as a tool for observers as well as clarify its ability to 
discriminate among models. Expectations for the shape of the density profile have been made which 
can explain the existence of a universal profile and support the NFW form (Evans & CoUett 1997; 
Padmanabhan, Cen, Ostriker, k, Summers 1996). Some of these depend on hierarchical clustering 
explicitly. 

The dependency of a universal density profile on hierarchical clustering is explored in this 
section. First, some of the possible forms of the density profile will be examined. After this, the 
results of fitting the density profiles to clusters formed hierarchically and non-hierarchically will be 
given. 



5.1. Introduction 

Over the years, a number of forms for the density profile of a system of collisionless particles 
has been suggested. Of particular relevance to the dark matter in galaxy clusters are the Hernquist 
profile and the Navarro, Frenk, and White (NFW) profile. Both have a limiting form in the near 
and far field with a smooth transition about some characteristic radius. 

The Hernquist profile (Hernquist 1990), suggested from observations of spherical galaxies, is 
given by, 

Pc r[rs + rY 

characterised by a length scale, r^, and an overdensity, 5. 

The work of Navarro et al. (1995, 1996, 1997) has built a case for a universal density profile of 
the similar form, 

Pc r{rs + rY 

Syer &; White (1998) (hereafter SW) claims that this form is a byproduct of hierarchical structure 
formation. It gives the generalised density profile form, 

p{r) rf 



' s 



(7) 



The exponents, a and (3, correspond to the exponential dependence of the density on radius in the 
near and far field. That is, for r <C r^, p oc and equivalently, for r ^ Vg, p on r~^. 
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The forms of Hernquist and NFW are particular cases of the SW density profile. For the 
Hernquist profile, a = I and (3 = 4. Correspondingly for the NFW profile, a = I and (3 = 3. It 
is the shape of the density profile in the inner radii which determines a. The inner radii, however, 
are affected strongly by the resolution of the simulation. Indeed, the cluster profiles generally 
span only two orders of magnitude in radius, the bulk in the outer halo. Any attempt to fit a 
smoothly varying curve over this span and then glean information about the near-field and far-field 
dependency of the profile is thus problematic. 

In order to circumvent the ambiguity in derived values of a and (3 due to the limited range 
of the fit, a non-smoothly varying form is required which is discontinuous in the first derivative at 
the 'knee' separating the nominal near and far fields. This suggests the following functional form 
should be fit, 

P{r) ^ f (5ar"", r < r^; 
Pc I Spr'/^, r>rs; 
Since the profile must be continuous at r^, 

Sp = (9) 

' s 

Simulations of hierarchically formed clusters agree that the NFW profile is an adequate de- 
scription of the dark matter density distribution. Cole & Lacey (1996) concludes this using a series 
of simulations of = 1 universes in which the initial density perturbation power spectrum is var- 
ied. Modelling clusters in a variety of cosmologies at very high resolution (256^ particles) Thomas 
et al. (1998) found the density profiles of the clusters to follow the form of the NFW profile. These 
simulations were N-body only. 

Similar to the situation being tested here, Huss, Jain, & Steinmetz (1998) describes a set of 
pure N-body simulations of an isolated cluster in which the amount of substructure is controlled 
via manipulation of the velocity dispersion of the particles. From this data set, they conclude that 
the NFW form maintains its universality in non-hierarchical scenarios. 

There is not as much agreement as to the details of the density profile. There are several 
predictions for the shape of the density cusp. They usually involve discussions of the growth of 
structure from a self-similar initial density distribution. An exception to this is the result of Evans 
&; Collett (1997). Looking at the stability of clusters to perturbations due to binary encounters, 
it shows that a = is a stable solution of both the Fokker-Planck and collisionless Boltzmann 
equations. 

Hoffman &: Shaham (1985) assumes a spherically symmetric halo collapse from an initial density 
perturbation with a power-law exponent of n to derive a radial dependency for the dark matter 
density of the form 

3(k+3) 

p{r) DC r "(^^ (10) 



- 15- 

for the inner cusp and p{r) cx in the outer halo. A similar result is found using a more detailed 
analysis in Padmanabhan et al. (1996). Using the assumptions of self-similarity and stability in 
the form of the cluster in the non-linear regime (that is, the morphology of a virialised object only 
scales with time), Padmanabhan et al. (1996) derives a relation for the 2-point correlation function 
in the non-linear stage of evolution from an initial linear density perturbation with a spectral index 
of n, 

3(n+3) 

^(r) ocr '(^. (11) 

Since the 2-point correlation function is simply the excess probability of finding any particle at a 
distance r from a given particle, if most of the matter is already located in the high-density regions, 
it is approximately related to the density profile by 

^-i + C(0 

Pbg 

This predicts 

for our density profiles. The same result is derived in a slightly different form in Padmanabhan 
(1996). Linear theory for the growth of clusters from an initial density fluctuation power spectrum, 
P{k) DC fc" predicts the same expression for a for p cx r" in the case of r <C r^. 

Extending the previous arguments to hierarchical clustering, Syer &: White (1998) argues that 
the initial cusp form found in the first objects formed will be maintained despite mergers. Consider 
a cluster absorbing a smaller satellite. If the density cusp in the cluster is steeper than that of 
the satellite, then the satellite will be tidally destroyed, softening the cluster cusp. If the satellite 
cusp is steeper, it will survive tidal disruption and sink to the bottom of the cluster, steepening 
the cluster cusp. 

Based on the self-similar infall model, which excludes hierarchical clustering, Henriksen &; 
Widrow (1998) finds p cx for the inner dark matter density profile and in the outer limits. 

Observational support exists for the NFW profile. Carlberg et al. (1997) claims agreement of 
the NFW profile with the average mass-density profile measured from 15 clusters in the Canadian 
Network for Observational Cosmology (CNOC) cluster survey. It finds 0.2 < < 0.3 for the 
average profile of their sample. However, it notes that the fit to a Hernquist profile has smaller 
residual errors. Others disagree. Adami et al. (1998) finds that the King profile, first suggested for 
globular clusters (King 1962), which has a core, fits the majority of a sample of clusters extracted 
from the ENACS (ESO Nearby Abell Cluster Survey) better than the NFW form. It notes, however, 
that this is not the case for every cluster in the sample; some are consistent with an NFW profile. 



(12) 
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5.2. The density profiles of the simulated clusters 

Dark matter density profiles were found for all the clusters formed hierarchically. These profiles 
are scaled readily by i?200 for each cluster to agree to a common form (Fig. 5). It is clear that there 
are three regimes, nominally the near-field, the far-field, and the background, each with distinctly 
different power law dependencies. Only the near and far fields concern us, here. 

The four forms of density profiles, NFW, Hernquist, SW, and the form suggested here (Eq. 8), 
were fit to the scaled SPH-estimated densities of the particles for the unsmoothed hierarchical run 
(Fig. 6). The method of fitting minimises x^- To separate the mean background density from the 
fitting, only the span r < 2i?2oo was fit. The parameters for the fits are outlined in Table 4. 

All forms fit the data equally well. However, due to the limited range in radii over which to 
fit the profiles (typically only two orders of magnitude), the limiting cases of the near and far field 
are poorly discerned. In particular, the value for a is poorly constrained. It is apparent that the 
'knee' position, Vg, and the near-field dependency, a, are being used as free parameters to shape 
the smooth transition between the two regimes. Consequently, neither represent the transition 
radius and near-field dependency one would expect. The SW fit fails entirely to discern an inner 
slope. However, for the outer slope, the SW fit, as well as the form suggested here, agree with the 
expectation of NFW. 

For these simulations, n = —1 which gives a = 1.5 using the predictions of Eq. 13. This 
compares quite well with the results of the fit to the form suggested. It also agrees with the results 
of fitting to the form suggested by Syer and White to the non-hierarchical data, for which the 
near-field is properly fit since the characteristic radius is further out. However, there is no reason 
that Eq. 13 should hold for non-hierarchical growth. 

Other groups find similar results. For a Cold Dark Matter (CDM) initial power spectrum, 
Moore et al. (1997) and Moore et al. (1998) find a = 1.4 in a high-resolution (softening of 5kpc) 
collisionless simulation. More generally, Fukushige & Makino (1997) finds for the core cusp a = 1 to 
2, in high-resolution collisionless simulations of a CDM model using the GRAPE special-purpose 
N-body hardware. The authors reported a resolution of 1 kpc in these simulations. 

The results of the fits for all the models (Table 5) indicate that the fits to the NFW form do 
as well as the less constrained form suggested here, even in the non-hierarchical cases. The scale 
radius increases with increased smoothing for the NFW fit in the high resolution simulation, but 
not the low resolution one. The inner profiles, as parametrized by a, have no dependency on the 
degree of smoothing but do exhibit a great deal of scatter. The mean value for the high resolution 
runs is (a) = 1.8 it 0.2. The outer profile does tend to become less shallow as the smoothing 
increases. This is qualitatively consistent with Crone, Evrard, & Richstone (1994) who find the 
profile becomes more shallow as the initial power spectrum index, n, is increased in a set of pure 
N-body simulations of clusters. 
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6. MASS-TEMPERATURE SCALING LAW 

In a hierarchical clustering scenario, the process of merging 'recreates' the cluster morphology 
on a continual basis. Clusters at a given epoch are being formed from structures created at a variety 
of earlier ages which are, in turn, formed by the amalgamation of many other smaller structures. 
In a non-hierarchical clustering scenario, the matter in the clusters has actually collapsed for the 
first time only recently. In both cases, if hydrostatic equilibrium has been established and there 
exist density and temperature profiles common to all clusters, then it can be expected that there 
will be a common mass-tcmperature scaling law for the clusters in a given clustering scenario. 

This universality of the gas density and temperature profiles has not, however, been established. 
It will be shown that the profiles affect the normalisation, r, of the mass-temperature relationship. 
Presuming the expression, T = tM , the steeper the profiles are, the lower the normalisation, r, 
is. 

In what follows, we will derive the mass-temperature scaling relation, making our assumptions 
clear. The normalisation factor will be more precisely defined, including its dependence on the shape 
of the gas density and temperature profiles which will then be examined. The actual normalisation 
factor for the simulated clusters will be presented, followed by a comparison with that determined 
implicitly from the derived relation combined with the information of the profile shapes. 



6.1. Derivation of the mass-temperature scaling law 



The equation of hydrostatic equilibrium relates the pressure, P{r), with the the mass, M(r), 
internal to the radius, r, via 

GM{r) _ 1 dP 
p dr ' 



(14) 



where p is the gas density and G retains its normal use as the gravitational constant. Combined 
with the perfect gas law, 

^ ^pT 

P=-^, 15 
p-ruu 

which relates the pressure to the temperature, T, with k , fi, and niu being Boltzmann's constant, 
the mean molecular weight, and the atomic mass unit, respectively, Eq. 14 gives, 

GM{r) 



kT f r dp ^ r dT\ 



r pruu \p dr ' T dr J 

Using the definition of In, the following relation is found: 



GM{r) kT{r) 



UrUu 



dlnp d In T 

H 

d\nr dlnr 



(16) 
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It is shown in Sections 6.4 and 6.5 that both the gas density and temperature profiles can be 
approximated by power laws away from the 'knee' radius. Taking the exponents of these respective 
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profiles to be Up and ht gives, 



p{r) oc r"'',r(r) oc r"^ 



GM{r) _ _kT{r) 



[np + nr] 



(18) 



In general, Up and ut are negative. 



The properties of clusters are characterised by a radius, Rg, and mass, Mj, of some overdensity. 
This overdensity is typically taken to be around the value of the virial overdensity for a top-hat 
collapse, i.e. 5 = 178. However, the characteristic radius is not particularly sensitive to the choice 
of 6. Given this description, we have the definition 

1/3 

(19) 



which, when combined with Eq. 18, gives 

nl/3 



47r 5pc_ 



G 



in ^ 



M 



2/3 _ kT{R-s) 



[rip + tit] . 
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Thus the scaling relationship, T oc M^/^, is found. Recall that the critical density is given by, 

_ 3(/il00km s-^Mpc-^)2 

On the scales of interest, this relation may be written as. 
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1 keV — {up + Ut) 



10^5 Mq 



(22) 



Significantly, the coefficient of the mass-temperature relation is dependent on the forms of the 
gas density and temperature profiles through the exponents Up and nr- 



6.2. The temperature normalisation parameter: r 

Since the temperatures of the clusters are scaled by the mass of the cluster, it is convenient to 
introduce the parameter 

r(r) 



r(r) = i^^^. (23) 

f _M200_\ 3 
Mq ) 

Thus, the r profile for a cluster will have the same form as T(r) but will be scaled by the cluster 
mass. Hence, if T(r) oc r"'^, then r(r) oc r"'^ as well. 

For ease of comparison, two temperature parameters are defined. The first, 

T(r=0) 
r = IkeV 

'num — 2 ! V / 

( M200 \ 3 
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can have its value found from the distribution of temperature versus mass for the clusters. 
Given that Eq. 22 gives T{Rs) while 

'Tnum is related to T(r — 0), we can compare the two by 
including a factor for the scaling between r = Rg and r = Riso, the radial limit of the isothermal 
core. This gives for the second temperature parameter, 

_ 2.83/.V3^-V3 

which will allow us to compare the actual central temperatures with the central temperature to 
be expected for a cluster with known density and temperature radial profiles at R^ which have 
power-law dependencies of Hp and ht respectively as well as isothermal radii approximated by 
Riso- Both of these normalisation parameters may be measured directly and independently from 
the simulated data. Their equivalence would verify the derived M-T scaling law encompasses the 
relevant physics. 



6.3. Mass-temperature distribution 

The mass-temperature relation was derived directly from the simulated data using the mean 
temperature in the core of the clusters and the mass within the virial radius, -R200- Temperatures 
within 0.2i?2oo are essentially constant (see Sec. 6.5). This distance was used as the cutoff radius 
in determining the mean cluster temperature, T^. The distribution of temperature versus mass is 
illustrated in Fig. 7. 

Fitting with a free coefficient for the relation T oc M", the correlation between central cluster 
temperature, T^, and the mass of the cluster is found to obey the relationship, 

since n = 0.67 ± 0.04, invariant of the degree of hierarchical clustering. The values for n for each 
simulation are listed in Table 6. 

For Table 6, errors were estimated using a bootstrap approach. Fits were calculated repeatedly 
using cluster data from a sample of half of the clusters, chosen randomly. The errors were deter- 
mined from the variance of the coefficients of these fits. This method, of course, biases towards 
the more common small-mass clusters whose counterparts in the real world would be less readily 
observed. When the samples are restricted to the more massive clusters, the results do not change 
significantly. Fits to the data points assume a correlation among clusters and, as such, lead to 
smaller errors than are found assuming independence among clusters (cf. Eq. 27). 

Using a fixed value of n = Vsj the normalisation factor, Tnum-, was calculated for the clusters 
within each simulation (Table 7) . The data points were assumed to be independent for the sake of 
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the fit and error estimate. That is, the error is, 



a = std 



Tr/kcV 



(27) 



/ _M200_\ 



2/3 



There may be a weak (Icr) dependence of Tnum on the degree of smoothing in the initial conditions. 
The value of Tnum drops from 6.8 ± 1.6 in the unsmoothed case to 5 it 1 for the 7Mpc smoothing. 
This is consistent with the trend in Tq found for the fits to r-profile. 

The value for the normalisation coefficient, Tnum, is in agreement with the number of numerical 
simulations, albeit on the high side. For the simulations of six clusters formed in a cold dark matter 
(CDM) model described in Navarro et al. (1995), a normalisation of 5.3±0.7 was found. Simulations 
of a variety of cosmological models which varied in are described in Evrard et al. (1996). There, 
a factor value of 5.20 is reported for r. Note that for this last value, the empirical translation 
M200 ~ 1-2M500 has been assumed here. Prom the data presented here, this ratio is found to be 
1.17 lb 0.11 and a slight trend towards higher values for more massive clusters is noted. For the 
previous assumption, this trend is insignificant. For the ten clusters described in Eke, Navarro, & 
Frenk (1998), a normalisation of 6 it 1 fits their data after scaling from the overdcnsity of 100 used 
in their analysis to the overdensity of 200 used here. These clusters were simulated using initial 
conditions for a low density, fiat universe. Using an Eulerian code, Bryan Sz Norman (1998) found a 
lower value of 4.7 it 0.1 from a variety of CDM models. Balogh, Babul, &; Patton (1998) has shown 
that if the gas is preheated and allowed to collapse adiabatically into an isothermal potential well, 
then the mass-temperature scaling relationship overestimates the halo masses by up to an order of 
magnitude. This occurs for haloes the size of groups of galaxies or less (M < 10^^ Mq ), which is at 
the low end of the size of the haloes examined here. Otherwise, it recovers the relationship found 
here, with a coefficient of about 4.5. 

The weighted mean of the value for tau using Tnum and the values mentioned above is 5.1 ±0.1. 
For this, an error of 5 times the least significant digit was attributed if none was otherwise given. 

The effect of galactic winds on the state of clusters is addressed in Metzler &; Evrard (1997). 
It reports that the gas density profiles are made more shallow by the winds while the central 
temperatures are not. The gas density profile parameter, Up rises from —2.34 to —1.75 with the 
inclusion of winds. As for the mass-temperature relation, for their sample of 18 clusters, the 
normalisation factor changes slightly from 4.8 ±0.4 to 5.2 ±0.3 with the inclusion of galactic winds. 
That the factor rises as the density profile becomes more shallow is consistent with Eq. 25 in trend, 
if not magnitude, since a value of 4.8 would rise to 5.5 if ht were held constant at —0.5. Metzler & 
Evrard (1997) does not report a value for nx though it indicates that the winds lead to a steeper 
temperature profile. This would partially offset the increase in r due to the decrease in Up. 

Various authors have found observational evidence supporting this scaling law. Schindler 
(1996) and Tsai, Katz, & Bertschinger (1994) find T = 7.8 keV (MTot/lO^^ M© ), with the cosmo- 
logical factors, h, set to that used in these studies (0.65). Horner, Mushotzky, & Scharf (1999) 
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verifies the form of the mass-temperature scaling relation using a sample of 23 clusters extracted 
from catalogues of clusters for which temperatures have been determined using ASCA (Fukazawa 
1997; Markevitch 1998) and masses via the projected velocity dispersion of member galaxies (Gi- 
rardi, Giuricin, Mardirossian, Mezzetti, &: Boschin 1998). It reports a normalisation factor of 



In order to calculate a value for Tanaiy, the gas density profile is required. A mean profile is 
determined for the clusters within a simulation using the gas density versus radius distribution 
for the cluster gas particles. The radii arc scaled by -R200 under the assumption that the gas will 
primarily follow the dark matter whose density distribution has been seen to scale with the virial 
radius, -R200- 

The gas densities were recalculated using the SPH density estimator with the smoothing lengths 
set to enclose exactly Ngpu = 32 particles. The distance from the cluster centres, scaled by 
the overdensity radius -R200 for each of the clusters, was found for each particle. The sample of 
particles comprised those that are within 2i?200 but beyond the point interior to which there were 
300 particles. This interior limit removes those particles for which the SPH density estimates are 
erroneous due to the steep density gradient (see Appendix A) . Fig. 8 illustrates the sample for the 
clusters formed hierarchically. The cluster profiles are scaled remarkably well by the virial radius, 
-R20O) considering that this is a combination of gas particles from 100 clusters spanning a range in 
mass of 30 times. It was verified that this scaling by the virial radius holds for the non-hierarchical 
simulations, also. A density profile was fit using the discontinuous density profile form given by 
Eq. 8 introduced in Sec. 5.2. The results are summarised in Table 8. 

There is a weak trend in which increases with increased smoothing. A much more significant 
trend is for a to decrease with increased smoothing. That is to say the density profile becomes more 
shallow as smoothing is increased. The outer profile changes little with a power-law dependency 
on radius of ^ —2.7 which is slightly shallower than the density profile for the dark matter (see 



The standard form to which gas density profiles are fit is the "/?" -model (Cavaliere & Fusco- 
Femiano 1976) which describes the expected density profile of an isothermal cloud. This model is 
used frequently in the interpretation of x-ray observations of galaxy clusters (Fabricant & Goren- 
stein 1983; Jones & Forman 1984). It has the continuous form 
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6.4. Mean gas density profiles 



Sec. 5.2). 
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The significance of j3 in this fit is that it also represents the ratio of the galaxy (or more specifically. 



collisionless baryonic component) kinetic energy to gas thermal energy. That is, 



P = 




(29) 



k T/ixrUu 



where Ugai is the one dimensional velocity dispersion for the galaxies in the cluster. Though there 
is no collisionless baryonic component in these simulations and the gas is well thermalized, it may 
still be interesting to see if the inferred value of (3 varies among the models. When Eq. 28 is fit to 
the mean gas profiles (Table 9), it is found that (5 varies little among the models with a mean value 
of 0.82 ± 0.03 for the high resolution runs. This compares with the value of 0.76 it 0.06 found by 
Tsai et al. (1994) who analysed a single cluster simulated at high resolution in a CDM model. Eke 
et al. (1998) simulated clusters in a low-density universe and found a mean value of 0.74 it 0.15. 



Like the gas density profile, the mean temperature profile is required to calculate Tanaly Since 
the cluster temperature is scaled by the mass, it can be expected that the temperature profiles 
are also accordingly scaled. As such, the profile of the scaled-temperature parameter, r, defined in 
Eq. 23, should be used to compare temperature distributions. 

A data set for the temperature parameter, r, was created in a manner similar to that of the 
gas density data set in Sec. 6.4. For the hierarchical case, the mean scaled r profile is shown in 
Fig. 9. There is a great deal more scatter than for the density. Some of this is due to smaller 
clusters being satellites of larger clusters. The halo gas of these smaller clusters gets shocked to 
high temperatures which then are scaled by the mass of the smaller cluster. However, there is an 
approximately isothermal core that extends to 0.2i?200- The outer profile may be crudely fit with 
a power-law dependency on radius. Fits to the profiles with the discontinuous form given in Eq. 8 
which has a free power-law index for the inner regime, a, find values of a spanning —0.05 to 0.10. 
For this reason, an approximate r-profile is fit using a form similar to the discontinuous form but 
with an iso-thermal core: 



where To is the scaled temperature (in units of T/M^^^) in the centre of the cluster. The results of 
fitting this to the scaled temperature versus scaled radius distributions for each of the simulations 
are given in Table 10. 

The central temperatures are found to generally decrease with increased smoothing. The 
isothermal core radii increase slightly with increasing smoothing, in step with the trend found for 
rg of the gas-density profiles, but marginally interior to r^. The power-law dependency on radius, 
riT, shows a trend to steepen slightly with increased smoothing, ranging from -0.4 to -0.6. This 



6.5. Mean r profile 




(30) 
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value has a strong dependency on resolution, with ut ~ —0.8 for the lower resolution runs, as does 
the iso-thermal radius which increases in diameter from 0.1 to O-SVa^^j, as the resolution is halved. 
This is consistent with the two-body heating problem described in Sec. 2.4. It was determined, 
there, that two-body interactions may be able to heat the gas to the temperature of the core over 
radii of 0.1 to 0.2 Mpc for the high resolution run and 0.1 to 0.5 Mpc for the lower. As such, the 
scaled temperature profile may become even more shallow with a further increase in resolution. 



6.6. Comparison with the semi-analytic prediction 

Recall Eq. 25, which defines, 

^ 2.83fe^/3^-V3 

Tanaly- + ^i.o" i^^J 

We can compare Tanaly to that inferred from the numerical simulations using np = (3 from Table 8 
and riT from Table 10. Table 11 summarises the results with the values of Tanaly given in the second 
column. 

Within error, tlie values for t from tlie sca,ling Iciw fit, '^num-) 

the T-profile fit, Tq, and the 

semi-analytic expectation assuming hydrostatic equilibrium, Tanaly , agree with the exception of the 
values for the 7 Mpc top-hat model. However, the trend for r to decrease with increasing smoothing 
is not as obvious when Tanaly is considered. In Sections 6.4 and 6.5, it is seen that Hp is invariant 
and riT is only weakly variant among the runs. It is the increase in Rig^ with increased smoothing 
that leads to the decrease in Tanaly (recalling that ut is negative). The slight increase in jit with 
smoothing tends to offset this trend. Consequently, within the scatter among clusters, the value 
for r does not vary appreciably and could be considered a constant among cosmologies. 



-24- 



7. Conclusions 

Using a series of numerical simulations, the importance of hierarchical clustering to the state 
of the matter in cluster of galaxies was examined. The simulations contained a coUisional compo- 
nent, permitting an examination of the hot gaseous halo of the clusters along with the dominant 
collisionless dark matter component. The comparison was carried out using simulations in which 
structure formed hierarchically, from smaller initial structures, and those in which the formation 
of these earlier structures was suppressed by smoothing the initial density distribution. 

The clusters within each simulation were not examined individually, but instead were examined 
as a whole sample. As required, individual cluster parameters were scaled by either the radius or 

mass at an overdensity of 200. 

Wc find the overall morphology of the simulated clusters robust to their formation method and 
relatively independent of the hierarchical clustering nature of structure formation. In particular, 
the dark matter profile and the mass-temperature relationship do not vary substantially between 
clusters formed hierarchically and those formed otherwise. However, these features do vary sys- 
tematically in small ways which may be significant to both theory and observations. These will be 
discussed in the following sections. 

It is noted that the gaseous component of the simulated clusters for this work were shown to 
be in hydrostatic equilibrium to beyond their respective virial radii, as indicated by R2oo- At this 
radius, the ratio of the actual gas pressure to the gas pressure expected, from the integration of 
the equation of hydrostatic equilibrium, is within 20% of unity for the majority of the clusters. 

7.1. Universal density profile 

The density profile form of Navarro et al. is found to fit the mean dark matter density profiles 
of the clusters of all models studied, indicating that the universal profile is not a by-product of 
hierarchical clustering. This profile form, however, fits only marginally better than the Hernquist 
profile. Indeed, the discontinuous form suggested here generally fits as well, if not better, than the 
NFW form. However, it has more free parameters. 

The discontinuous form indicates that the inner regimes of the dark matter density profile 
are approximately dependent on radius with p oc r^^'^, independent of the degree of smoothing 
of the initial conditions. The outer regimes, in contrast, do seem to be dependent on smoothing. 
The profile becomes marginally shallower as the degree of smoothing is increased. The dependency 
varies from p oc r~^'^ in the case of the hierarchical, unsmoothed case to p oc r^^'^ for the 7 Mpc 
top-hat smoothed model and the (14Mpc)~^ frequency cutoff model. The characteristic radius, 
rg, is affected by the smoothing; it expands from 0.1i?2O0 to 0.2i?2oo when hierarchical clustering is 
eliminated. This is consistent with the later epoch of formation of the cluster. 
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The NFW profile has been established to be a robust description of the dark matter density 
profile, independent of hierarchical clustering and, by extension, independent of the assumed cos- 
mology. This strongly implies that dark matter in galaxy clusters is distributed in this manner. If 
detailed mass maps of clusters using, for example, gravitational lensing observations indicate that 
the mass is not distributed in this fashion, it would indicate that either numerical simulations are 
missing some significant physics or the Standard Model is flawed in a fundamental manner. It 
will be difficult to attribute the discrepancy to the assumed initial density perturbation spectrum 
within the Standard Model. Unfortunately, this also implies that the observed mass density profiles 
of clusters will tell us little about this initial perturbation spectrum. 



7.2. Mass-temperature scaling law 

The clusters were found to have iso-thermal cores. Using the temperatures within these cores 
and the masses within an overdensity radius of 200, a normalisation factor, Tnum, for each simula- 
tion was found. Integrating the equation of hydrostatic equilibrium permitted a derivation of the 
temperature parameter, Tanaly, which is dependent on the radial profiles of both the gas density 

and the temperature as well as the radius of the iso-thermal core. Mean profiles of these parame- 
ters were determined for the clusters within each simulation. This permitted the determination of 
values for Tanaly which were compared with the respective values of Tnum to both verify the model 
and explain any variation found between the models. 

While determining gas density profiles for calculating Tanaly i the value of (3 was also determined. 
This factor is related to the ratio of galaxy kinetic energy to baryonic thermal energy in real systems 
as well as being the /? used in fits of the gas distribution in the data from observations of clusters. 
It is found to not vary with the degree of smoothing. Instead, it maintains a value of 0.82 it 0.03. 

The mass-temperature scaling law (Eq. 26) is dependent on the degree of hierarchical cluster 
formation, albeit at a la level. The temperature parameter Tnum varies from 7itl in the unsmoothed 
case to 5±1 in case of the 7Mpc top-hat run with the values for the intermediate smoothing models 
lying in between. The trend is consistent with the expectation of hydrostatic equilibrium. This value 
decreases from 7.7 it 0.2 to 6.4 it 0.7 over the span of smoothing. There is a systematic discrepancy 
between the expectation, Tanaly, and the value derived directly from the clusters, Tnum, with the 
semi-analytic expectation consistently larger than the measured value. This may be attributed to 
numerical, two-body heating which is found to be significant in the domain of Riso- But again, the 
discrepancy is less than the chistcr-to-clustcr variation within a simulation. The results of the fits 
to the profiles indicates that it is the variation in the isothermal radius that explains the variation 
with smoothing. The isothermal radius increases as the degree of hierarchical clustering is reduced. 
This occurs in step with the gas density profiles which become shallower in the inner radii. 
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A. SPH NEAR STEEP DENSITY GRADIENTS 

It is known that the SPH estimation of the density for a particle fails near a steep density 
gradient. The density is essentially determined by finding the radius of a sphere that encompasses 
some given number of particles, Nsph, which is typically 20-30. That the particle contributions 
arc weighted by a smoothly-varying kernel function is a second-order correction to the estimate. 
Crudely, in the presence of a varying density gradient, contributions from particles in the higher- 
density region exceed those from the lower-density region, leading to an over-estimate of the local 
density. Non-spherically symmetric kernels alleviate this problem, but are computationally expen- 
sive both in memory resources and CPU time. In this appendix, this phenomenon will be explored 
in more detail with the aim to predict those gradients that suffer this effect most strongly. 

For any spherically symmetric weighting kernel, the SPH estimate of the density for a particle 
in a constant density gradient will produce no error. This requires, of course, that the gradient be 
constant over the span of x — 2h < x < x + 2h which implies that the gradient is resolved by the 
particle density. Here, we will examine the error induced by a varying density gradient. 

The toy model that will aid us in our exploration will be a volume with an increasing density 
gradient in one direction. Specifically, the density will behave as p oc x~^. This is a form similar to 
the density profiles found in cosmological objects. An iso-density surface for these objects, though 
spherical in morphology, is locally flat. Hence this toy model is relevant to the density gradients 
found in cosmological objects. 

The models were constructed from 100 slices, each with the same number of particles but with 
widths, Wi, such that Wi-iF = Wi = Wi+i/F where F >1. To test variation with form, F is set 
here to 1.1 and 1.5 for two different sets of gradient fields. The particles were distributed randomly. 

In our tests, the mass per particle is unity, giving a number density, n, equivalent to the mass 
density, p. 

Error in the density estimate arises from both statistical variance as well as the systematic error 
due to the density gradient. The former error is reduced by using more particles in the estimate of 
the density. The latter is reduced by kernel averaging over a smaller region, and hence using fewer 
particles. The regime of dominance of these effects was explored by calculating the density using 
^SPH = 64, 32, 16, and 8 neighbour particles. Near the density cusp, where gradients are higher, 
the density estimate which averages over fewer particles (Nsph ^ 16) displays large variance, but 
no systematic over-estimation of the density. In comparison, the estimates using larger numbers 
of particles {Nsph > 16) systematically over-estimate the density by a factor of 10 to 20% within 
X = 0.1 of the cusp but exhibit Vs less variance. Further away from the cusp, where the density 
gradient is smaller, the larger values of Nsph allow density estimates with smaller amounts of 
statistical variance and no systematic bias. 

Increasing the number density of particles increases the resolution of the simulation. This 
permits steeper gradients to be resolved. For the p oc x~^ profile of our tests, this behaviour was 
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qualitatively verified by varying the number of particles per slice. Numbers of 50, 100, 200, 500, 
and 1000 particles per slice were used. For the gradient regime, > 500, the smaller number 
densities increased the error in the density estimate. The gradient field that used 1000 particles 
per slice was able to go to an order of magnitude greater density gradient. 

More useful is a parametrization of the density gradient that is independent of the local number 
density of particles. Consider the parametrization based on the change of the density over the 
smoothing distance, 2h, given by p{x — h)/p{x + h). It was found that this parametrization is 
independent of the number of particles per slice. It is also independent of the form of the gradient. 
It indicates that SPH fails (at the 20% level) when the density across a smoothing length varies by 
more than a factor, Fc, or 3-4 times, for the p oc profile used here. 

An estimate of the radial distance from the centre of a dense clump at which SPH will begin 
to calculate erroneous density values can be made using Fc. This factor, F, can be defined as 

„ n(r — h) , . , , 

F = -) -( (Al) 

n(r + h) ^ ^ 

where n is the number density of particles, or p/m. The critical radius, Tc, below which the density 
estimates will be in error is related as 

n{rc — h{rc)) 
nivc + hire)) = 

The estimate n(r + Sr) n{r) + Sr^ gives 

n{r - h) ^ n{r)-h{r)^ 

n{r + h) n{r) + h{r)^' ^ ' 

If the density profile is assumed to take the form, n(r) = nor~°' then dn/dr = —anor~°'~^ leading 
to 

n(r — h) nor~°' + ah{r)nor " ^ 1 + ah{r)r ^ r + ah{r) 



n(r + h) UoV — ah{r)nor " ^ 1 — ah{r)r ^ r — ah{r) 
Using the relation between h, Nsph, and n, 

1 /3 iVgpj, y/^ 



(A4) 



and Eq. A2, the value for rc can be derived, 

i-f^a/3iV^y/^^ 

2\A TTUo J Fc-1 ^ ' 

In these studies, Nsph = 32 and, as seen, Fc ~ 3, giving the approximate equation, 

1-^ 2a 

rc ' ^ -Tj^- (A7) 
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This allows the derivation of an approximate minimum radius for a given cluster with a compression 
factor (c.f. overdensity) of Sc at the critical radius, Vc- The compression factor is defined as the 
density at the critical radius compared with the mean density in the volume. That is, 

- (A8) 



where the notation (v4)y refers to the volume- weighted mean of the parameter A. Since (n) 
i?es'^, where Res is the number of nodes per side in the initial density distribution (i.e., A', 



V 



gas 



Res^), 

ric = ScRes^. (A9) 



Recall that n(r) = implying Ug = ScRes^r^ giving, 

2a 



SJ^Res 



(AlO) 



For the maximum overdensity observed, 5c = 10^, Tc = 40kpc for simulations with 64^ particles. 
For 5c = 10^, this increases to Tc = 400 kpc . Clearly, rc is enclosing a constant number of particles. 
This minimum number of particle, A^mm = ^Trr^ScRes^ , when rc from Eq. AlO is substituted, gives, 

N^i^c^4{2af. (All) 

This number is the minimum number of gas particles, within rc, required to properly resolve the 
density gradient given by a. For a = 2, the required number is 250 but this number is sensitive to 
a and drops to 100 for a = 1.5 

The number, Nmin, was used to set the minimum cluster size in the analysis as well as limit the 
regime over which bulk properties were determined; radii enclosing fewer than N^in gas particles 
were ignored. 
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Ho 65 km s-^Mpc-^ 

^DarkMatter 0-9 

^Gas 0.1 

A 

Power law index —1 

as .935 



Table 1: Properties of the assumed cosmology. 





^ of particles 


rUg 


rriDM 


'^smooth 


e 




gas 


dark 


IQio 


Mo 


(/i-^Mpc) 


(/i-^kpc) 


unsmoothed 


643 


643 


1.03 


9.27 





20 


3Mpc top-hat 


643 


643 


1.03 


9.27 


3 


20 


(7Mpc)-^ k-cutoff 


643 


643 


1.03 


9.27 


7 


20 


7Mpc top-hat 


643 


643 


1.03 


9.27 


7 


20 


(14Mpc)-i k-cutofF 


643 


643 


1.03 


9.27 


14 


20 


low-res. unsmoothed 


323 


323 


8.24 


74.16 





40 


low-res. 7Mpc top-hat 


323 


323 


8.24 


74.16 


7 


40 



Table 2: Properties of the simulations. Given are the number of each type of particle, the mass per 
particle for both the gas, rUg, and dark matter ,mDM, the effective smoothing radius, r smooth , and 
the gravitational softening length, e. The length of a side of the simulation volume is 40h~^ Mpc 
for all simulations. 





■^clusters 






unsmoothed 


100 


0.05 


1.52 


3 Mpc top-hat 


48 


0.05 


2.85 


(7Mpc)-i k-cutoff 


48 


0.05 


1.75 


7 Mpc top-hat 


8 


0.03 


1.95 


(14Mpc)-i k-cutoff 


6 


0.15 


2.27 


low-res. unsmoothed 


3 


0.54 


0.70 


low-res. 7 Mpc top-hat 


5 


0.40 


1.88 



Table 3: Results of the cluster search. Given for each run is the number of clusters found, Ndusters-, 
and the range of masses of the clusters, Mmin and M^ax- 
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rs a (3 x^/N 

NFW 0.06 ±0.01 1 3 0.98 

Hernquist 0.19 ±0.01 1 4 1.00 

SW 0.0008 ±0.003 -7 ±3 2.94 ± 0.02 0.98 

here 0.12 ±0.01 1.84 ±0.03 2.74 ± 0.01 0.97 



Table 4: Coefficients for fits to the dark matter density profiles for the hierarchical case using a 
variety of profile forms. The SW fit, though free to fit an inner slope, fails to converge to a value 
for a and fits, instead, a single power law. 



Rim 


f'\nTll 




(1 


3 






NFW 


0.06 ± 0.01 


1 


3 


0.98 


unsmoothed 


here 










0.12 ±0.01 


1.84 ± 0.03 


2.74 ± 0.01 


0.97 


3Mpc top-hat 


NFW 


0.10 ±0.01 


1 


3 


0.98 






here 


0.07 ±0.01 


1.5 ±0.2 


2.46 ± 0.01 


0.98 


(7Mpc)-i 


NFW 


0.09 ± 0.01 


1 


3 


0.97 


k-cutofF 


here 


0.3 ±0.2 


2.1 ±0.2 


2.8 ±0.2 


0.99 


7Mpc top-hat 


NFW 


0.28 ± 0.01 


1 


3 


0.99 






here 


0.23 ±0.01 


1.70 ± 0.03 


2.23 ± 0.01 


0.98 


(14Mpc)-i 


NFW 


0.21 ±0.01 


1 


3 


0.99 


k-cutoff 


here 


0.18 ±0.01 


1.78 ± 0.03 


2.28 ±0.01 


0.98 


low-res. 


NFW 


0.26 ± 0.03 


1 


3 


1.01 


unsmoothed 


here 


0.23 ± 0.01 


2.08 ± 0.06 


2.38 ± 0.05 


1.01 


low-res. 


NFW 


0.24 ± 0.02 


1 


3 


1.00 


7Mpc top-hat 


here 


0.26 ±0.1 


4±2 


2.4 ±0.1 


1.40 



Table 5: Coefficients for fits to the dark matter density profiles for all the models. 



'Tnum 



unsmoothed 


5.57 


± 


0.50 


0.61 


±0.03 


3Mpc top-hat 


5.89 


± 


0.61 


0.67 


±0.04 


(7Mpc)~^ freq. cut 


6.96 


± 


0.63 


0.75 


±0.04 


7Mpc top-hat 


4.80 


± 


0.86 


0.67 


±0.12 


(14Mpc)~^ freq. cut 


6.96 


± 


0.63 


0.75 


±0.04 


low-res. unsmoothed 


5.7 


± 


0.5 


0.66 


±0.07 


low-res 7Mpc top-hat 


5.7 


± 


1.5 


0.68 


±0.26 



Table 6: Mass-temperature scaling law with free coefficients. 
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unsmoothed 


6.8 ± 1.6 


3Mpc top-hat 


6.2 ± 1.1 


(7Mpc)-^ k-cutoff 


6.2 ± 1.1 


7Mpc top-hat 


4.9 ±0.9 


(14Mpc)-^ k-cutoff 


5.5 ±1.1 


low-res. unsmoothed 


5.8 ±1.3 


low-res. 7Mpc top-hat 


4.9 ± 1.0 



Table 7: Mass-temperature scaling law coefficient. The coefficient, t„„„j, of Eq. 26 is given for the 
simulations. It relates the cluster temperature (in keV ) with the cluster mass (in 10^^ M0 ). 













a 






(3 






unsmoothed 


0.43 


± 


0.03 


2.37 


± 


0.01 


2.45 


± 


0.01 


0.96 


3Mpc top-hat 


0.56 


± 


0.03 


2.27 


± 


0.01 


2.80 


± 


0.04 


0.93 


(7Mpc)-i k-cutoff 


0.54 


± 


0.01 


2.24 


± 


0.01 


2.97 


± 


0.02 


0.94 


7Mpc top-hat 


0.59 


± 


0.01 


1.60 


± 


0.02 


2.99 


± 


0.03 


1.00 


(14Mpc)-i k-cutoff 


0.23 


± 


0.01 


1.78 


± 


0.03 


2.34 


± 


0.01 


0.94 


low-res. unsmoothed 


0.43 


± 


0.05 


3 


± 


2 


2.34 


± 


0.03 


1.00 


low-res. 7Mpc top-hat 


0.73 


± 


0.03 


1.54 


± 


0.08 


3.0 


± 


0.1 


1.01 



Table 8: Coefficients for fits to the gas density profiles. All use the discontinuous form. 







P 




unsmoothed 


0.020 ±0.001 


0.808 ± 0.001 


0.96 


3Mpc top-hat 


0.044 ± 0.001 


0.840 ± 0.002 


0.94 


(7Mpc)-i k-cutoff 


0.056 ± 0.001 


0.868 ± 0.002 


0.95 


7Mpc top-hat 


0.145 ±0.005 


0.811 ± 0.001 


1.04 


(14Mpc)-i k-cutoff 


0.073 ± 0.002 


0.791 ± 0.004 


0.94 


low-res. unsmoothed 


0.02 ± 0.02 


0.78 ± 0.02 


4.05 


low-res. 7Mpc top-hat 


0.42 ± 0.04 


1.05 ± 0.05 


1.02 



Table 9: Coefficients for the /3-fit to the gas density profiles. 
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To 




—riT 




unsmoothed 


8.5 ±1 


0.10 ±0.03 


0.41 ± 0.04 


1.02 


3Mpc top-hat 


8.6 ±0.1 


0.10 ±0.03 


0.47 ±0.01 


1.01 


(7Mpc)-i k-cutoff 


7.52 ± 0.04 


0.16 ±0.01 


0.57 ±0.01 


1.03 


7Mpc top-hat 


4.85 ± 0.04 


0.26 ±0.01 


0.56 ± 0.02 


1.01 


(14Mpc)-i k-cutoff 


6.59 ± 0.04 


0.18 ±0.01 


0.57 ±0.01 


1.06 


low-res. unsmoothed 


9.9 ±0.5 


0.5 ±0.1 


0.7 ±0.2 


1.03 


low-res. 7Mpc top-hat 


4.0 ±0.1 


0.58 ±0.02 


0.92 ± 0.07 


1.00 



Table 10: Coefficients for fits to the scaled gas temperature profiles. All use a discontinuous form 
with an iso-thermal core. 





Tanaly 


Tnum 


To 


unsmoothed 


7.7 ±0.2 


6.8 ± 1.6 


8.5 ±1 


3Mpc top-hat 


6.8 ±0.2 


6.2 ± 1.1 


8.6 ±0.1 


(7Mpc)-i k-cutoff 


7.1 ±0.3 


6.2 ± 1.1 


7.52 ± 0.04 


7Mpc top-hat 


6.4 ±0.7 


4.9 ± 0.9 


4.85 ± 0.04 


(14Mpc)-i k-cutoff 


7.1 ±0.3 


5.5 ±1.1 


6.59 ± 0.04 


low-res. unsmoothed 


7.6 ±0.4 


5.8 ±1.3 


9.9 ±0.5 


low-res. 7Mpc top-hat 


5.4 ±0.4 


4.9 ± 1.0 


4.0 ±0.1 



Table 11: Analytic and numeric coefficients for the M-T relation. Given are the semi-analytic 
estimates for the scaled temperature parameter, Tanaly^ the numerical value from fitting to the 
temperature-mass distribution, Tnum, and the values derived from fitting the mean profiles, Tq. 



-36- 



10 



-4 



10" 



0) 

o 

CL 



10 



-6 



io; 



1 1 1 

^^^^^ ^ 


1 1 


II II 

unsmoothed 
3 Mpc top-hat 

- - (7 Mpc)"'' k-cut 
7 Mpc top-hat 

__. (14Mpc)"^ k-cut 


1 






1 


\ ^^^^ 
\ \ 

\ 

^\ 




1 1 1 


K _ 

' ■ \ 
» \ 
I 1 \ 
r. 1 \ 

!■. , \ 

1 ■ ■ \ 

1 \ 

1 •. 1 \ 

1 ■ 1 \ 

1 ■ 1 I I 1 1 



Mpc 



-1 



0.02 0.03 0.04 0.06 0.08 0.1 0.15 0.2 0.3 0.4 



Fig. 1. — The initial power spectra of the simulations. 
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Fig. 2. — The projected gas densi- 
ties of the simulation boxes at t = 1 
{i.e., present). Top left plate is evolved 
from the unsmoothed initial condition 
(IC); middle left from the 3Mpc top- 
hat smoothed IC; middle right from the 
(7Mpc)-^ cut-off filtered IC; bottom 
left from the 7Mpc top-hat smoothed 
IC; and bottom right from (14Mpc)'~^ 
cut-off filtered IC. Darker shades indi- 
cate higher densities, in log scale. 
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Fig. 3. - Comparison of the gas pressure profile in the simulated clusters to that expected for a 
hydrostatic scenario for those clusters formed hierarchically. Plotted are the pressure profiles for 
the clusters with the three highest central pressures (left) and the three lowest central pressures 
(right). The solid line is calculated from the integration of the equation of hydrostatic equilibrium 
while the dotted line is the measured gas pressure profile. The vertical line denotes the virial radius, 

-^200- 
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Fig. 4. — Comparison of the measured gas pressure profile to that expected for a hydrostatic 
cluster formed non-hierarchically. Plotted are the pressure profiles for six of the seven clusters in 
the sample. As in Fig. 3, the data for the three clusters with the highest central pressures are on 
the left, and those for the three lowest central pressures are on the right. The solid line is calculated 
from the integration of the equation of hydrostatic equilibrium while the dotted line is the actual 
gas pressure profile. The vertical line denotes -R20O) the approximate virial radius. 




Fig. 5. — Dark matter density profiles for the clusters formed hierarchically at the epoch corre- 
sponding to the present. The radii have been scaled to R2qo- 



10^ 




Fig. 6. — Hierarchical data. Fits to the mean dark matter density profile for the forms of Navarro, 
Frenk, and White (top left), Hernquist (top right), Syer and White (lower left), and the form 
suggested here (lower right). The dashed lines are the fits. The solid lines are the mean profiles 
with la variations bound by the dotted lines. 
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Fig. 7. — The temperature versus mass for the clusters. The sohd hne is the fit T = 5.9M^/^. 
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Fig. 8. — Gas densities for the hierarchical run. The densities are for the sample of particles 
described in the text. The radii have been scaled to i2200- The total sample corresponds to about 
30% of all gas particles in the simulation. Contour lines are of constant mass surface density in 
r/R2Q0 — T space, separated by factors of 10 Z^. 
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Fig. 9. — The scalcd-tcmperature parameter, r, for the hierarchical run. The sample of particles 
is the same as in Fig. 8. Contour lines are of constant mass surface density in r/i?200 ~ t space, 
separated by factors of 10 . 



